Ultrahigh energy cosmic rays as heavy nuclei from cluster accretion shocks 
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Large-scale accretion shocks around massive clusters of galaxies, generically expected in the cold 
dark matter scenario of cosmological structure formation, are shown to be plausible sources of the 
observed ultrahigh energy cosmic rays (UHECRs) by accelerating a mixture of heavy nuclei including 
the iron group elements. Current observations can be explained if the source composition at injection 
for the heavier nuclei is somewhat enhanced from simple expectations for the accreting gas. The 
proposed picture should be clearly testable by current and upcoming facilities in the near future 
through characteristic features in the UHECR spectrum, composition and anisotropy, in particular 
the rapid increase of the average mass composition with energy from 10 19 to 10 20 eV. 
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Introduction. Several decades after their discovery, the 
origin of UHECRs with energies 10 18 -10 20 eV and above 
remains one of the biggest mysteries in physics and astro- 
physics [H, 0] . The observed global isotropy in the arrival 
directions strongly suggests that they are of extragalactic 
origin. However, no unambiguous identification with any 
type of astrophysical source has been achieved so far. 

If UHECRs are protons, photopion interactions with 
the cosmic microwave background (CMB) should induce 
severe energy losses at > 7 x 10 19 eV for propagation 
lengths > 100 Mpc [3|. It is also possible that UHE- 
CRs are dominated by heavy nuclei, particularly at the 
highest energies where the mass composition is observa- 
tionally quite uncertain Q. In this case, photodisinte- 
gration and photopair interactions with the far-infrared 
background (FIRB) and the CMB govern the energy loss 
length Q, which is ~ 300 Mpc at 10 20 eV for iron nuclei, 
somewhat larger than that for protons Q . 

Only a few types of astrophysical objects appear capa- 
ble of accelerating UHECRs up to the highest observed 
energies, such as the jets of radio-loud active galactic 
nuclei (AGNs) and gamma-ray bursts (GRBs) [2J. Clus- 
ter accretion shocks have also been proposed as UHECR 
sources Q- In the currently favored hierarchical scenar- 
ios of cosmological structure formation, massive clusters 
of galaxies should be surrounded by highly supersonic 
accretion flows that give rise to powerful and long-lived 
shock waves extending over Mpc scales. Numerical simu- 
lations show that such shocks should possess high Mach 
numbers implying that they can accelerate particles 
with hard spectra at high efficiency through the first or- 
der Fermi mechanism [9(. However, the maximum energy 
attainable for protons seem to fall short of 10 20 eV by 1-2 
orders of magnitude 2, 11. 12|. 



sition of heavy nuclei, in particular Fe nuclei up to en- 
ergies ~ 10 20 eV. Note that "accretion shocks" here sig- 
nify not only shocks due to infall of cold gas from the 
diffuse intergalactic medium, but also those associated 
with gas inflow along large-scale filaments that are impor- 
tant in terms of energy dissipation 8] . Our findings hold 
for plausible assumptions on the source properties and a 
wide range of extragalactic magnetic field strengths, and 
should be clearly testable in the near future. 

Model. We adopt the cosmological parameters h = 0.7, 
O m = 0.3, £1a = 0.7 and og = 0.9. First we estimate 
the maximum energy by considering a fiducial, Coma- 
like cluster of total mass M = 2 x 10 15 M m , and extend- 
ing the discussion of Rcf. 



rjj to include heavy nuclei. 



In this letter we show that cluster accretion shocks may 
provide a viable explanation of the observed properties 
of UHECRs if they accelerate a suitably mixed compo- 



The accretion shock radius is taken to be the virial ra- 
dius, R s ~ 3.2 Mpc, so that the shock velocity relative to 
the upstream gas is V s = (4/3)(GAf/i? s ) 1 / 2 ~ 2200 km/s 
We choose B s = 1/xG for the magnetic field at the 
shock, as suggested by some recent observations around 
the virial radius [lj|. The timescale for shock accelera- 
tion of particles with energy E and charge Z is i acc = 
2Qk{E)JV? = {20/3)(Ec/ZeB s V s 2 ), assuming a parallel 
shock [42j, together with the Bohm limit for the diffu- 
sion coefficient k(E) which can be induced by CR wave- 
excitation [lH and is compatible with supernova remnant 
(SNR) observations [161 ] . The maximum energy -E ma x can 
be estimated by comparing i acc with timescales for limit- 
ing processes, such as energy loss by photopair and pho- 
topion interactions with the CMB for protons, and addi- 
tionally photodisintegration interactions with the FIRB 
and CMB for nuclei, as defined in Ref. [f|. Another limit 
may be the diffusive escape time from the acceleration re- 
gion, t esc ~ R 2 /5k(E) [ljj, which we can equate with i acc 
to obtain E nl&x /Z ~ (3eB s R s V s /10c) ~ 6.5 x 10 18 eV. Fi- 
nally, the Hubble time ta = 1.4 X 10 10 yr represents an 
absolute upper limit. These comparisons are illustrated 
for selected species in FigQ] For protons, E max < 10 19 



2 



eV, limited either by photopair losses or escape, confirm- 
ing previous findings. However, heavy nuclei with higher 
Z and correspondingly shorter i acc can be accelerated 
to higher energies, with Fe reaching ~ 10 20 eV in the 
same conditions, notwithstanding photonuclear interac- 
tion losses. 




10" 10 Z ' 
log E [eV] 

FIG. 1: Comparison of timescales versus particle energy E 
at cluster accretion shocks for shock acceleration t acc (di- 
agnonal lines), and for energy losses from interactions with 
background radiation fields (curves), for protons (thick dot- 
ted), He (thin dotted), O (thin solid) and Fe nuclei (thick 
solid). Also indicated are the Hubble time tn (dashed line) 
and the escape- limited maximum energies (circles). 



For clusters of mass M, the rate of gas kinetic en- 
ergy flow through accretion shocks can be estimated as 
L acc ~ 9 x 10 45 (M/10 15 M o ) 5 / 3 erg s" 1 H[T|. We take 
the number density of clusters with M > 10 15 M Q to be 
n s = 2 x 10 _6 Mpc -3 , in accord with that observed locally 
within 120 Mpc 18], and consistent with theoretical mass 
functions to within a factor of 2 [h|. The power density 
is P acc - L acc n s = 2 x 10 40 erg s^Mpc" 3 . Around 10 19 
eV, protons should dominate with propagation lifetime 
~ 4 x 10 9 yr. To account for the observed energy density 
of UHECRs > 10 19 eV of - l(r 19 erg cm" 3 , the required 
power density is then P>i9 ~ 3 x 10 37 erg s _1 Mpc~ 3 , ne- 
glecting the effect of magnetic fields. Thus cluster accre- 
tion shocks should reasonably accommodate the energet- 
ics of UHECRs. 

For a quantitative study of this scenario, we simulate 
trajectories of nuclei above 10 19 eV, accounting for all rel- 
evant energy losses and deflections by extragalactic mag- 
netic fields (EGMF), including secondary nuclei arising 
from photodisintegration 2^, 21 1. Since the true nature 
of EGMF is currently very uncertain 22j, we consider 
different cases that may bracket the range of possibili- 
ties, specifically, models of relatively str ong EGMF that 



are not included, however. The observer is chosen in a 
low EGMF region with some resemblance to the Local 
Group's actual environment. Assuming n 8 and P acc as 
above, we simulate a large number of realizations where 
the locations of discrete sources follow the baryon den- 
sity, and obtain the average and cosmic variance of the 



spectra and composition [21 , 

A fraction /cr of the accretion luminosity L acc is 
converted to UHECRs with energy distributions oc 
E~ a exp(— E/E max ), with E max and normalization be- 
ing different for each species. We set E max /Z = 5 x 10 18 
eV, which is a fair approximation to estimates obtained 
as in FigfTJ The escape of UHECRs into intergalactic 
space is assumed to be energy-independent, at least in 
the limited energy range 10 19 -3 x 10 20 eV. This could be 
the case for diffusive escape in directions away from the 
filaments, or possibly advective escape during merging 
events. 

For the elemental composition at injection, we follow 
Ref. [1^1 in taking the abundance ratio by number of He 
to protons to be 0.042. All heavier elements are assumed 
to have the same relative abundances at fixed energy per 
nucleon as that inferred for the sources of Galactic cos- 

and scaled with respect 



mic rays at GeV energies [27|, [25 
to protons by the metallicity of the accreting gas ( rela- 
tive to the solar abundance. We take £=0.2 as suggested 
by both observations [5§| and theoretical simulations [1(3] 
for the warm gas falling in from filaments. This assumes 
that particle injection from the thermal gas at cluster ac- 
cretion shocks works in a way similar to Galactic cosmic 
ray sources, presumably SNRs, which is plausible consid- 
ering that the shock velocities are of the same order, and 
the temperature of the accreting preshock gas at 10 5 -10 6 
K is that of the hot phase of the interstellar medium. 
On the other hand, the nonlinear modification of shock 



structure by CR pressure in the upstream region [3l| can 
be stronger in cluster shocks than in SNRs owing to the 
acceleration to much higher energies, and this may pos- 
sibly lead to a further enhancement of heavier nuclei at 
injection [32j]. We account for this by an additional fac- 
tor A@ in the injected abundance of nuclei with atomic 
number A. In such conditions, progressively heavier nu- 
clei reach higher -E max and become the dominant species 
in the source spectrum, analogous to that inferred for 
Galactic cosmic ray sources in the knee region 

Results and Discussion. Fig[2] displays our results for 
the observed spectrum and composition for the case of 
a=1.7 and /3=0.5. The spectrum is consistent with the 
current data for the HiRes 33J and the Pierre Auger ex- 
periments [34| within about 2 sigma of combined cosmic 
and statistical variance, for both the differential and in- 
tegral fluxes. The latest analysis of the AG AS A data 



trace large-scale structure as in Refs. [2l|, [23[, as well 
as the case of no EGMF, more in line with some other 
models 24\. The effect of Galactic magnetic fields [25} 



35j] (not plotted) may possibly result in similar spectra 
depending on the choice of air shower simulations [3rjj ]. 
Note that values of a < 2 are naturally expected at the 
high energy end for nonlinear shock acceleration that ac- 
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counts for the dynamical back reaction from CRs 31 1. 
For (3, acceptable fits to the spectrum are obtained in 
the range 0.3-0.5 corresponding to a factor of ^3-7 en- 
hancement for Fe nuclei. It remains to be seen whether 
this can be explained from first principles in the context 
of nonlinear acceleration theories [HJ . We predict a clear 
steepening in the spectra > 10 20 eV, due not only to pho- 
tonuclear losses during propagation, but also because our 
estimated E max at the source is not far from this value. 
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FIG. 2: Observed UHECR spectrum (top) and mean mass 
composition (bottom) versus energy E (1 EeV = 10 18 eV) 
from cluster accretion shocks for the case of a = 1.7 and 
/3 = 0.5, compared with the current data for HiRes (bars) 
and Auger (asterisks). The histograms are the average result 
over different model realizations for the cases with (thick) 
and without (thin) extragalactic magnetic fields, and the thin 
curves outline the median deviations due to cosmic variance 
for the former case only. The straight line in the top panel 
denotes the injection spectrum. 



Normalizing to the observed flux, the required source 
power above 10 19 eV is inferred to be in the range 
P>iq ~ (1 - 50) x 10 38 erg s _1 Mpc~ 3 for the cases with 



EGMF, fluctuating by a large factor depending on the 
specific realization of the source locations with different 
propagation lengths to the observer. Comparing with 
Pace estimated above, we deduce a reasonable range of 
fcR — 0.005 — 0.3. We note that for a < 2 as expected in 
nonlinear shock acceleration, the energetics is dominated 
by CRs at the high energy end, while those at lower ener- 
gies do not contribute much power. For the cases without 



EGMF, P>i 9 ~ 3xl0 37 erg s^Mpc -3 , quite independent 
of the realization, leading to /or — 0.002. Low values of 
/cr could reflect inefficient escape from the system rather 
than inefficient acceleration. 

Consistent with claims from HiRes [37| , the flux below 
~ 2 x 10 19 eV is dominated by lighter nuclei. This is due 
both to their prevalence at lower energies in the source 
spectrum as described above, and to pileup of photodis- 
integrated fragments coming down from heavier nuclei at 
higher energies during propagation. The latter effect is 
more prominent with EGMF because of the longer prop- 
agation lengths. The rapid increase in the average mass 
composition above ~ 2 x 10 19 eV, a natural consequence 
of higher Z nuclei extending to higher E max , is a defini- 
tive prediction of our scenario to be verified by the new 
generation experiments. 

Since massive clusters are relatively rare in the local 
universe within ^300 Mpc, the UHECR flux should be 
dominated by only a few nearby sources. Nevertheless, 
strong deflections of the highly charged nuclei in EGMF 
allow consistency with the observed global isotropy of 
UHECRs with the current statistics. On the other hand, 
we predict that with a sufficient number of accumulated 
events, significant anisotropy should appear toward a 
small number of individual sources. Fig[3] compares the 
angular power spectrum of arrival directions expected for 
100 and 1000 events above 4 x 10 19 eV, in the case of 
a particular realization with a single, dominant cluster 
at distance D ~ 50 Mpc, somewhat resembling Perseus 
or Coma. The sky distribution may show an excess of 
events extended by few tens of degrees, offset from the 
true position of the cluster by a similar amount. Such 
trends may in principle reveal important information on 
the EGMF. Note, however, that these inferences could 
be affected when Galactic magnetic fields are included. 
For the case with no EGMF, the anisotropy is too strong 
to be compatible with the current data, although includ- 
ing weak EGMF and/or Galactic fields may alleviate this 



25j,|38|. At any rate, a general feature of scenarios with 
heavy nuclei as UHECRs is that a small number of pow- 
erful sources can be consistent with the current global 
isotropy and lack of association with known objects. This 
is in contrast to AGN scenarios with protons that must 
rely on sources more numerous than powerful FRII radio 
galaxies, since the latter are as rare as massive clusters 
despite being the most promising AGN candidates [17] . 

A good match to the spectra can also be obtained with- 
out the A$ factor if £ ~0. 5-1.0, but such high metallicities 
may not be expected for the accreting gas Si. The 
results are fairly insensitive to variations in a alone, since 
prescribing the source injection composition at fixed E/A 
implies a weighting factor A 01 ^ 1 at a given E [26| . which 
counteracts the E~ a dependence. On the other hand, a 
modest increase of E maK /Z by a factor of 2 produces a 
large spectral bump at a few times 10 19 eV from photo- 
disintegration products, at odds with the observations. 
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FIG. 3: Angular power spectrum C(l) of UHECR arrival di- 
rections above 4 x 10 19 eV versus multipole I, for a realiza- 
tion with extragalactic magnetic fields and a single, dominant 
cluster at D ~ 50 Mpc. The crosses are for 100 events with 
AGASA + SUGAR exposure and diamonds for 1000 events 
with Auger North + South exposure. Vertical bars indicate 
statistical errors. 



Thus it is crucial that the source spectrum cuts off near 
the photodisintegration threshold on the FIRB, which is 
naturally expected in our cluster scenario (FigQ}, but is 
not obvious for nuclei acceleration in other contexts, e.g. 
AGNs [13 or GRBs H- 

Although most previous studies of UHECR nuclei con- 
centrated on the propagation aspects 28 1, here we 
have presented a physically plausible scenario of UHECR 
nuclei production by cluster accretion shocks that can 
account for the current observations. Combined with 
complementary information from X-ray and gamma-ray 



observations [12|, |39| , detailed measurements of UHECR 
composition and anisotropy with facilities such as the 
Pierre Auger Observatory, the Telescope Array, and the 
future Extreme Universe Space Observatory should pro- 
vide a clear test of whether the largest bound structures 
in the universe are also the largest and most powerful 
particle accelerators. 

After this paper was submitted, the Auger collabora- 
tion announced new observational results regarding the 
composition [i(3] and anisotropy 41 1 of UHECRs. Our 
prediction of heavy-dominance at the highest energies 
is concordant with the composition estimates, whereas 
it could be disfavored by the anisotropy results. Al- 
though the implications of these apparently conflicting 
results are not yet clear, it is certainly possible that some 
fraction of UHECRs are heavy nuclei, which can be ac- 
counted for by the scenario discussed in this paper. 
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